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Quantum state tomography JJ, the ability to deduce the state of a quantum system from measured 
data, is the gold standard for verification and benchmarking of quantum devices. It has been realized 
in systems with few components [2HZ]: but for larger systems it becomes infeasible because the 
number of quantum measurements and the amount of computation required to process them grows 
exponentially in the system size. Here we show that we can do exponentially better than direct state 
tomography for a wide range of quantum states, in particular those that are well approximated by 
a matrix product state ansatz. We present two schemes for tomography in 1-D quantum systems 
and touch on generalizations. One scheme requires unitary operations on a constant number of 
subsystems, while the other requires only local measurements together with more elaborate post- 
processing. Both schemes rely only on a linear number of experimental operations and classical 
postprocessing that is polynomial in the system size. A further strength of the methods is that the 
accuracy of the reconstructed states can be rigorously certified without any a priori assumptions. 



I. INTRODUCTION 

One of the principal features distinguishing classical 
from quantum many-body systems is that quantum sys- 
tems require exponentially many parameters in the sys- 
tem size to fully specify the state, compared to only lin- 
early many for classical systems. Put to use construc- 
tively, the exponential complexity enables the construc- 
tion of information processing devices fundamentally su- 
perior to any classical device. At the same time, how- 
ever, this "curse of dimensionality" makes engineering 
tasks — such as verifying that the quantum processing de- 
vice functions as intended — a daunting challenge. 

The full determination of the quantum state of a sys- 
tem, known as quantum state tomography has al- 
ready been achieved by measuring a complete set of ob- 
servables whose expectation values determine the quan- 
tum state PHZ]- As it is typically formulated [5], simply 
to output an estimate for a generic state would take ex- 
ponential time in the system size iV, given that there are 
an exponential number of coefficients in a generic state's 
description. This is but one of several inefficiencies. Most 
quantum states have exponentially small amplitudes in 
almost every basis, so to distinguish any one of those am- 
plitudes from zero, one must take an exponential number 



of samples. Assuming one were able to collect all the data 
from an informationally complete measurement, one is 
left with the intractable computational task of inverting 
the measured frequencies to find an estimate of the state. 

However, the traditional representation of quantum 
states is in a sense too general. Indeed, states which oc- 
cur in many practical situations are specified by a small 
number of parameters. An efficient description could be 
a practical preparation scheme which outputs the state; 
or, in the case of thermodynamical equilibrium states, a 
local Hamiltonian and a temperature. This insight is not 
new: researchers in many-body physics and quantum in- 
formation theory have found many classes of states which 
are described by a number of parameters scaling polyno- 
mially in N [5141 3| and which closely approximate the 
kind of states found in physical systems [HI [THj . How- 
ever, the question of whether these restricted classes can 
be put to use in the context of tomography has remained 
largely open. 

In this work, we address the above challenge. The 
physical system we have in mind is one where the con- 
stituents are arranged in a one-dimensional configura- 
tion (e.g., ions in a linear trap But the methods 
that we are presenting here can be generalized to higher- 
dimensional arrangements such as those realized in opti- 
cal lattices TF]. It is highly plausible that in such a set- 



ting, correlations between neighboring particles are much 
more pronounced (due to direct interaction) than corre- 
lations between distant systems (mediated e.g. by global 
fluctuations of control fields). An efficiently describable 
class of states anticipating exactly this behavior has long 
been studied under the names of finitely correlated states 
(FCS) or matrix product states (MPS) 0110]. Impor- 
tantly, restricting to this class is not a limitation since 
every state may be written as a MPS with a suitable, al- 
beit possibly large, matrix dimension. Since many states 
that are relevant for quantum information processing or 
quantum many-body physics have a small (independent 
of N) bond dimension, our methods are directly applica- 
ble to such states; examples include, but are not limited 
to, ground and thermal states, the GHZ, W, cluster, and 
AKLT states, the latter two being universal resources 
states for quantum computing. 

Given that standard tomography is no longer feasible 
in a range of recent and upcoming experiments involving 
large numbers of qubits, our results represent a signifi- 
cant advance in the ability to verify and quantitatively 
and efficiently benchmark systems of experimental im- 
portance. 



II. RESULTS 

In the following we present two schemes for identify- 
ing systems that are well approximated by an MPS, ini- 
tially focusing on pure states for simplicity. We will view 
each system as consisting of a linear chain of N qudits, 
each having dimension d. Both schemes require the mea- 
surement of linearly (in the system size N) many local 
observables within finite accuracy, polynomial classical 
post-processing of the data and can certify the accuracy 
of the reconstructed state without making any technical 
assumptions about the state in the laboratory. The first 
scheme requires unitary control and local measurements 
while the second scheme removes the need for unitary 
control at the cost of more elaborate post-processing. 



A. Scheme based on unitary transformations 

The key idea of the method consists of a sequential 
procedure to disentangle the left hand side of the chain 
from the right hand side, using a sequence of unitary 
operations with small interaction length independent of 
A'^. The result will be a product state and a sequence 
of local unitary operations from which to construct the 
original state. 

Suppose the ideal state in the laboratory is g = 
which we assume for clarity is a MPS of given bond di- 
mension. This implies that the rank of reductions to 
one part of a bipartite (left vs. right) split of the chain 
is bounded by a constant R. The protocol starts by 
estimating, through standard state tomography, the re- 
duced density matrix of the first k = [log^(i?)] -I- 1 sites. 



\4>) 



J L 



T" 

|o) 



Ui 

1 r 



U2 



T 

|o) 



1 r 



1 r 

|0) : 




FIG. 1. Quantum circuit that transforms \(l>) into 
|Qji®]v-K+i ^ ^j^j^ K = 3. The unitaries tji successively 
disentangle the particles and the state |r;) on the last sites 
acts as a boundary condition to determine the MPS descrip- 
tion of |(^). 



a « trK+i....,Ar[p]. This reduced density matrix has the 

eigendecomposition a = X]^=i '^r\4>r){4>r\ where the rank 
R < d"^^. Hence there exists a density matrix with one 
fewer qudit that has the same rank R and eigenvalues cr^ 
as a. Therefore we can disentangle the first site in a with 
the following unitary acting on the first k sites: 



d-l < 

E 

s=0 



s'=0 



(1) 



where . . . , \4'r) have been extended in some arbi- 
trary way to get a complete basis for sites 1, . . . , k. Ap- 
plying U produces the state U\<j)) = \^)i®\v)2,...,n , where 
\v) is some pure state on sites 2, . . . ,N. Hence U disen- 
tangles the first qudit from all the others. Now, set aside 
this first qudit, look at sites 2, . . . , k + 1, and repeat the 
above process as indicated on Fig. [T] In this way, one ob- 
tains a sequence of unitaries Ui, . . . ,U, 



N-K+l, 



where each 



acts on sites i, . . . ,i + k — 1. This sequence transforms 
10) into Un-.+i ■■■W) = |o)®JV-«+i ^ 1^^^ ^j^g^g 
is some pure state on the last k—1 sites. 

In summary, this scheme infers the quantum circuit 
used to prepare an MPS llTj . The MPS decomposition of 
10) can then be obtained readily from the Ui and \r]) [18]. 
If the state in the laboratory g is arbitrary, then the 
reduced density matrices a will generally have full rank. 
Hence in each step we will need to truncate the k qudit 
state (7 to a rank R subspace with R < d'^~^. Then the 
above method will produce an MPS approximation to g. 
The accuracy of this estimate can be certified, without 
any assumptions on the state, by keeping track of the 
effects of truncating each of the reduced states a. As 
shown in the Methods, errors of magnitude e due to finite 
measurement precision or truncation error (as measured 
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by the weight of the truncated space) accumulate at most 
linearly with the number of sites, and can be evaluated 
directly from the data, leading to an estimate of accuracy 
Ne. 

The present scheme requires unitary control of k neigh- 
boring qudits, which is challenging to implement in many 
current experimental settings. We now present a second 
scheme that avoids unitary control and requires only lo- 
cal measurements. 



B. Scheme based on local measurements: 
Certification of estimated state 

Consider again the state g in the laboratory and sup- 
pose that a tomographically-complete set of local mea- 
surements on all groups of k neighbouring qubits has 
been performed. Further suppose that sufficient data has 
been taken to yield estimates ai of the local reductions 
Qi = t'Ci^,„^i.^+k+i,...,N[g\ such that 

lift - CTilltr < Q • (2) 

Determining these approximate reduced density opera- 
tors ffi completes the experimental work. In the remain- 
der we describe the classical post-processing that will re- 
sult in a MPS estimate 1-0) to g and a lower bound to the 
fidelity {iplgl-i/j) |19| that does noi require assumptions on 
the nature of the state g. 

We start with the latter, since the following easy 
calculation yields the fidelity bound we have in mind 
and hints at the MPS estimate we are after. Suppose, 
for the sake of the argument, that is the unique 
ground state (with energy zero) of a local Hamiltonian 
H = ^^hi, where the hi is a projection operator act- 
ing only on sites i -\- 1, . . . ,i -\- k (as it turns out, generic 
MPS are of this type). Then, expanding in the eigenbasis 

2" — ! 

tv[Hg] > AEY,{En\g\En) = ^E{1 - mW), (3) 

n>0 

where we denoted by AE the energy gap above the 
ground state 1-0). Hence, we have the fidelity bound 

(Vl^l^) > 1 - ^^IM > 1 - ^ ( E tr[k-^] + ^0 • 

i 

(4) 

In other words, the Hamiltonian acts as a witness for its 
ground state \ip). This bound is tight: suppose the ex- 
perimental estimates ai and the reductions of the state 
in the laboratory gi coincide, that is Ci = 0. If, in addi- 
tion, the reductions of li/j) match the ai then, as was 
assumed to be the unique ground state with energy zero 

oij^t^i, we have X^i t^'I^^'^'i] = 0; i-S-, {tp\g\i^) = 1- 

The goal is now clear: find a local gapped Hamilto- 
nian H such that the reductions of its ground state are 



close to the ai. A priori it is unclear whether such a con- 
venient witness always exists and how it could be found 
efficiently. However, using formal methods, one can show 
that if the true state g is close to a generic MPS, then 
such a witness Hamiltonian always exists [10 . What is 
more, it can be constructed from the estimate of the algo- 
rithm sketched below. Its properties, chief among them 
the gap, are efficiently computable. In the Methods sec- 
tion, we detail the efficient computation of these quan- 
tities, and we also consider how to treat states such as 
the GHZ for which local marginals alone are not quite 
sufficient for complete characterization (they violate our 
"generic" condition). 



C. An illustrative example 

We illustrate how our certification procedure operates 
if our estimate for the state p is a linear cluster state 
PU] . The cluster state is defined as the unique eigenstate 
(with eigenvalue -1-1) of stabilizers Ki — Zi^iXiZi+i 
for all i = 2, . . . , — 1 (together with boundary terms 
X1Z2 and Zn^iXn, which we do not treat separately for 
simplicity). Assume that we have performed standard 
quantum state tomography on sets of three neighbour- 
ing spins, k — 3, which is the smallest useful set because 
in a cluster state the rank of the reduced density ma- 
trices of contiguous blocks is upper bounded by i? = 4. 
Let us now assume that on the basis of these data, our 
procedure suggests that the linear cluster state is indeed 
our estimate. The local Hamiltonian in this case is given 
by H = X;,(l - Ki)/2, where the (1 - K^)/2 = h, are 
projectors, H has the cluster state as its unique ground 
state (with energy zero) and an energy gap AE = 1. The 
fidelity between the cluster state \ipcs) and the state g is 

bounded by {ipcs\9\ipcs) = 1 - Y.^{^Af^'^^i\ + where 
Ci from Eq. (|3| quantifies the statistical and experimen- 
tal error in the local experimental estimates, and tr[ft,i(7j] 
quantifies how much the laboratory state g deviates from 
an exact cluster state. 

Tomography on a cluster state can also be performed 
with the scheme based on unitary transformations. A 
cluster state is the output of a quantum circuit where 
each qubit is initially prepared in the state |-|-) = 
l/\/2 (|0) + |1)) and a controlled phase transformation 
CZ acts successfully on each pair of neighboring qubits. 
The CZ gate changes the sign of state |11) and acts 
trivially on the other states of the computational basis. 
Thus, a cluster state is the output of a circuit whose 
structure corresponds to the one indicated on Fig. [T] with 
K = 2 and the scheme based on unitary transformations 
directly applies. Note that the unitary scheme takes ad- 
vantage of the decreased rank of a reduced density matrix 
on the boundary to save one qubit worth of local tomo- 
graphic effort (k = 2 vs. = 3). 
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D. Scheme based on local measurements: Efficient 
determination of an MPS estimate 

With the experimentally obtained ct^, we now turn to 
the task of finding an MPS 1-0) such that its reductions 
tri,...,i;i+fe+i,...,7v[|'0)(''/'|] closely match the di. In other 
words: Let be all possible products of Pauh oper- 
ators (enumerated by m) that act non-trivially only on 
sites z + l,...,i + /c. Then, as the are an orthogonal 
basis, the cTi may be expanded as 



2N 



(5) 



where the expectations tr[(TP^] G R are obtained as re- 
sults of tomographic measurements. Then we need to 
find a matrix such that for all m and i the ex- 

pectations tr[|'(/;)('0|P^] coincide with those of the tomo- 
graphic estimates, tr[|'(/') (V'|Pm] — tr[CTiPm]. 

The method of choice for such a problem is singular 
value thresholding (SVT) [3TJ |^ , which has been devel- 
oped very recently in the context of classical compres- 
sive sampling or matrix completion |23) and may also be 
applied to the quantum setting [24W28]. SVT is a re- 
cursive algorithm that provably converges to a low-rank 
matrix satisfying constraints of the type tr[|-0)(-0|P^] = 
tr[(TiP^]. Unfortunately, a straightforward application of 
SVT requires time and memory that scale exponentially 
with the number of particles. However, a modification of 
the algorithm allows us to overcome this problem. 

Given the measured values tr[(TiP^] the algorithm may 
then be described as follows. First set up the operator 
R = ■ tr[(7iP4]P^/2^ and initiahze Yq to some arbi- 
trary matrix (e.g., the zero matrix). Then proceed recur- 
sively by finding the eigenstate |j/„) with largest eigen- 
value, Um of Yn and set 



Yn+1 —Yn + Sn{R — X„). 



(6) 



So far, this algorithm still suffers from the fact that in ev- 
ery step the 2^ x 2^ matrix Yn needs to be diagonalized. 
However, the Yn are of the form J2mi '^m,iPm: '^rn,i G R, 
where the P^ act non-trivially only on sites . . . , i+k, 
i.e., they have the form of a local "Hamiltonian" . Hence, 
\yn) can be determined as the highest energy state of 
this Hamiltonian. For this task standard methods have 
been developed in condensed matter physics [131 for 
which the number of parameters scale polynomially in 
the system size and converge rapidly [30l |3T] to the op- 
timal MPS approximation. The standard but exponen- 
tially inefficient SVT algorithm possesses a convergence 
proof while our efficient modification does not. Hence, 
we now present numerical examples for different target 
states 10) to demonstrate the feasibility and efficiency of 
the proposed algorithm. 
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FIG. 2. Performance of the MPS-SVT algorithm for the 
ground state (the target state of the critical Ising model. 
We chose k = 2, i.e., only nearest neighbour reductions are 
used to reconstruct the state. Plot illustrates the scaling of 
the error in the fidelity 1 — /]v,n = 1 — K0l3/n)P ~ N"^ /n with 
the number of spins and the iterations n of the algorithm. 



These numerical examples suggest convergence of our 
algorithm to a MPS that closely matches the experimen- 
tally obtained reductions ai. To arrive at the fidelity 
bound we follow the steps of Section II B (i) Obtain es- 
timates of the reductions to k adjacent spins Qi of the 
state in the laboratory such that ||(7i — ^i||tr < e^, (ii) 
compute the expectations p^. = tr[(7iP,'„], which are the 
input to the MPS-SVT algorithm, (iii) obtain a MPS es- 
timate \yn)^ the reductions of which closely match the 
di by utilizing the MPS-SVT algorithm. As \yn) is an 
MPS, one can then efficiently obtain a parent Hamilto- 
nian [TU] and a lower bound. A, on the energy gap above 
the ground state (see Methods for details). Putting all 
this together, the above programme returns a state |y„), 
its parent Hamiltonian H — hi, and a number A such 
that 



{yn\Q\yn) > 1 - 



Y^ii^i + tr[/ii(Tj 



(7) 



E. Example: Strongly interacting quantum systems 

We start with ground states of nearest-neighbor Hamil- 
tonians on a chain, i.e., the |0) = \gs) are completely de- 
termined by all the reductions to two adjacent spins and 
we expect the above algorithm not only to produce states 
that match the reduced density matrices of the ground 
states but, in fact, states that are themselves close to the 
ground states. Among ground states of one-dimensional 
nearest-neighbor Hamiltonians, those at a critical point 
are the most challenging to approximate by MPS as they 
violate the entanglement area law ^32J. We demonstrate 
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FIG. 3. Fidelity /]v,n = K0|j/ii)p5 for a- fixed value n = 5, as 
a function of the number of spins A'^ for the ground state (the 
target state 10)) of the random nearest-neighbour Hamiltoni- 
ans of Eq. [S] As in Fig. [2j we chose fc = 2 to reconstruct the 
state. The plot shows densities (arrows indicate the mean) 
obtained from 1000 random realizations. Similar to the Ising 
model, the scaling for fixed n of 1 — fN,n is better than . 

the effectiveness of our algorithm for such an example: 
the critical Ising model in a transverse field on a chain 
of length N with open boundary conditions. The ground 
state of this Hamiltonian is unique and hence it is com- 
pletely characterized by its reductions to fc = 2 neigh- 
bouring sites. In other words, if we find an MPS that has 
the same reductions, the fidelity will be one. We proceed 
as follows, (i) We obtain the reductions Qi to two neigh- 
bouring sites of the true ground state, (ii) from these 
reductions we obtain the expectations p„i. = tr[pii-'mj, 
which are the input to the MPS-SVT algorithm. In Fig. [5] 
we show the fidelity f^.n = |(yn|0)P of the true ground 
state \(p) and the n'th iterate of the above algorithm as 
a function of n and the length N of the chain. It shows 
that for fixed system size 1 — fN,n decreases as ^ 1/n 
while for fixed n it increases slower than iV^. This pro- 
vides an indication that our algorithm is polynomial in 
the system size. 

The Ising model is solvable and, in order to show that 
we are not considering a special case that is particu- 
larly favourable, we also consider one-dimensional ran- 
dom Hamiltonians of the form 

i? = E^N;'i, (8) 

where the , ■rl_^^ act on spin i and i + 1, respectively, 
and are hermitian matrices with entries that have real 
and imaginary part picked from a uniform distribution 
over [—1,1]. For each Hamiltonian, as before, we first 



determine the ground state \gs) (our target state \(j>)) 
and its reductions and then computed the fidelity fN.n 
after n iterations of the MPS-SVT algorithm, see Fig. [3] 

F. Example: W-state preparation in ion traps 

Our method is of interest for many situations in which 
standard tomography will not be feasible. This is the 
case for the verification of state preparation in experi- 
ments with too many particles. An example is the recent 
ion trap experiment [S] for the preparation of W-states, 
10) = (|10---0)-H|010---0)-h--- + |0---01))/A/]V, that 
were limited to 8 qubits principally because the classi- 
cal postprocessing of data became prohibitive for longer 
chains. Here we demonstrate the efficiency of our ap- 
proach (we are not limited to few ions and demonstrate 
convergence for up to 20 ions - even higher number of 
ions are easily accessible due to the MPS alteration of 
the SVT method) by illustrating how one would postpro- 
cess experimentally obtained reduced density matrices to 
guarantee the generation of |0) or a state very close to it. 
Note that, as in the above spin chain examples, one need 
only take tomographic data on pairs (fc = 2) of neigh- 
bouring qubits. We mimic experimental noise by adding 
Gaussian distributed random numbers with zero mean 
to the Prrii ■ After initializing the algorithm with Yq = R, 
where we obtain R from the MPS representation of 10), 

weusea;„ := X^m, \P7n,i- {yn\P^^\yn) \ as a figure of merit 
for convergence, i.e., after a given number of iterations, 
we pick the \yn) with minimum a;„. The result of such a 
procedure is shown in Fig. [4] 

III. DISCUSSION 

We have presented two schemes that efficiently produce 
an MPS description as a tomographic estimate of a quan- 
tum state, along with a tight fidelity bound. We empha- 
size that no assumptions are necessary for our scheme; if 
no such MPS description exists, this will be evident from 
the local tomographic data and our schemes will herald 
a failure. However, the enormous successes of MPS for 
describing both one-dimensional quantum systems found 
in nature as well as a host of states relevant to quantum 
information ensures that our methods will be very useful 
in practice. 

So far we presented the method for pure states and one- 
dimensional systems. Various generalizations are pos- 
sible: Our first scheme using unitary control can also 
handle mixed states corresponding to small ensembles 
of pure MPS. Suppose we are presented with a state g 
that is a mixture of M pure states, each of which is an 
MPS having bond dimension D. Then the reduction of 
g to any subsystem has rank at most MD. We can pro- 
ceed as before, performing unitary operations on blocks 
of K = [log^(MD)] -I- 1 sites, in order to disentangle in- 
dividual sites from the rest of the chain. At the end of 
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FIG. 4. Performance of the MPS-SVT algorithm for W- 
states, 10) = (|10---0) + |010---0) + --- + |0---01))/yiV. We 
are not limited to few ions and demonstrate convergence for 
up to 20 ions - even higher number of ions are easily acces- 
sible due to the MPS alteration of the SVT method, demon- 
strating the efficiency of our approach. We mimic experimen- 
tal noise by adding Gaussian distributed random numbers 
with zero mean to the local expectations p^. — tr[(Ti_PmJ and 
show results for n = 4000 MPS-SVT iterations. Plot shows 
/]v,n = K0|j/n)P as a function of the number of ions, A'^, for no 
noise (dots) and Gaussian noise (densities obtained from 100 
realizations for each A'', arrows indicate mean) with a stan- 
dard deviation of 0.005 (blue; only even values of A'^ plotted 
for clarity) and 0.01 (black, odd A''). 



the chain, we will find a mixed state rj of rank M on the 
last K — I sites. We decompose as a mixture of M pure 
states. This yields a representation of g as an ensemble of 
M pure MPS, each with bond dimension at most dMD. 

Our second scheme using local measurements can also 
be extended to handle mixed states. While the /c-particle 
reduced density matrices do not uniquely determine the 
mixed state g, reconstructions of better and better qual- 
ity can be obtained by increasing k. As an example, 
suppose g is the Gibbs state corresponding to a /c-local 
Hamiltonian H , i.e. the state minimizing the free energy 



iT:[gH]-TS{g). 



(9) 



The first term is, as before, determined by the reduced 
density matrices. The entropy of the total state how- 
ever can only be learned exactly from the complete den- 
sity matrix. Yet, for essentially all reasonable physical 
systems, the entropy density limfc_j.oo '5'(trfc+i_...(^))//c in 
the thermal state of a Hamiltonian exists ^33i (In partic- 
ular, finitely correlated states are precisely those states 
whose entropy density is exactly equal to S{tT:k+i,...{Q))~ 
»S'(trfc,...(^)) for some finite value of k). As a consequence, 
the total entropy of the state can be estimated efficiently 
from knowledge of the reduced density matrices. Hence, 
our second scheme using local measurements may be ex- 



tended to mixed states by considering purifications and 
can then, e.g., also handle thermal states of local Hamil- 
tonians, under the physically reasonable assumption that 
the entropy density exists [33] . In addition, it can be gen- 
eralized to all mixed finitely correlated states (FCS) [S], 
though it is not always possible to certify the resulting 
estimates. 

Higher-dimensional systems are more challenging, be- 
cause the most straightforward generalization of MPS, 
known as projected entangled-pair states (PEPS) |34j . 
cannot be computed as efficiently. However, the certifi- 
cation method using a frustration-free parent Hamilto- 
nian remains efficient in the case of qubits with nearest- 
neighbor couplings |35j . Combinations of our techniques 
can be used to reconstruct other classes of states, such 
as tree tensor networks [36 and multiscale entanglement 
renormalization ansatz (MERA) states [37] , for which ef- 
ficient heuristics for minimizing local Hamiltonians are 
available. 
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IV. METHODS 

We start by recalling the MPS representation of a state 
I?/)) with open boundary conditions (generalisations to pe- 
riodic boundary conditions are entirely straightforward). 



di-l 

E- 

si=0 



Mi[si]---Mn[sn]\si 



■ sn) 



(10) 



siv=0 



where the Mjs] are Di x Di+i matrices with Di 
Dn+1 — 1- We denote the bond dimension by D 
max Di . 



A. Direct tomography 

This method proceeds by disentangling all the qudits 
of the chain sequentially. Thus, it yields a valid MPS 
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description if every unitary exactly disentangles one qu- 
dit. Put another way, while it is crucial to obtain a good 
estimate of the d'*" ^-dimensional subspace on which a is 
supported, it is not necessary to identify the eigenvectors 
of (7 exactly : any set of orthonormal vectors generating 
the subspace is sufficient for our tomography procedure 
and leads to an MPS description in another gauge [TU] , 
This property will be central to our error analysis. 

To understand the effect of errors and imperfections in 
our tomography procedure, consider the very first step 
of the recursive procedure. Tomography is performed on 
the first k sites to ideally find a state with non-maximal 
support, and unitary Ui is applied to rotate that state 
into the subspace W^^°^ = |0) ® (C^)®""!. In any ex- 
perimental setting, the resulting state tJi\4i) would surely 
not lie entirely in that subspace. This can be either be- 
cause the state of the system is not exactly an MPS with 
bond dimension D, but merely close to one or because 
our estimate of the density matrix a on sites 1 to k is 
slightly wrong due to measurements that are, in prac- 
tice, noisy and restricted to finite precision. Indeed, we 
can expect that the reduced density matrix on the first k, 
sites will actually be full-rank, though most of its prob- 
ability mass will lie on a subspace of dimension at most 
D. So, each time we apply a disentangling operation t/^, 
we also want to truncate the reduced state to the desired 
subspace. Similarly, a faulty estimate of a will result in 
a small probability mass that lies outside the estimated 
support. 

Given an estimated disentangling unitary C/i, any state 
\4>) can be expressed as 



be 



|0)^l^i) + |ei) 
Vl + (ei|ei) 



(11) 



where |ei) is some sub-normalized error vector supported 
on the subspace orthogonal to The unitary tJi 

is chosen to minimize the norm ei = (ei|ei) of this error 
vector. It is possible to directly estimate the magnitude 
of the error ei by measuring the first qudit in the stan- 
dard basis; the error is equal to the population of the 
non-zero states. 

In subsequent steps of the recursion, we are given a 
state of the form 



U,...UA(h) = 



|0)^' ® W) + 

v/1 + (erier 



(12) 



where |e™^) is the accumulated error vector that lies in 
the subspace orthogonal to = |0)®* ® (^^d^®N-z^ 

As a first step, we can truncate this error vector by mea- 
suring the first i particles in the standard basis and post- 
select on the all-zero outcome. This occurs with a proba- 
bility roughly 1 — e"", and leaves the system in the state 
|0)®* (g) \rii). We then repeat the steps leading to Eq. 



11 



with the post-selected state \rji). The resulting state will 



C/,+1 [/.... C/i 



|0)®*^C/,+i|7y,) + C/,+i|e? 



VI + (erkD 



|0>' 



\m+i) 



where 



/„cm I cm \ 



+ |er) 



(13) 
(14) 

(15) 



and therefore ef^^ < -I- e^+i. 

Thus, we see that errors accumulate linearly with 
the number of particles; if we denote = 
Ul... ?7j^_„+i|0)^-'"+i|r/Ar_^+i) the estimated MPS, we 
have 



N-K + 1 

< E (16) 

i=l 



where e — max^ |||ei)||. The overall error is at most the 
sum of the individual errors on each step. In addition, 
each of the is revealed during the tomographic proce- 
dure because they correspond to the post-selection suc- 
cess probability. This provides a direct method to certify 
the inferred state. 



B. Parent Hamiltonians 



Let be as in ( 10 1 and such that J2s ^'^i i'^W^ i^V = ^ 
for all i — 1, . . . , A^. This can always be achieved by 
subsuming qudits at the beginning and end of the chain 
into qudits with higher dimension and successive singular 
value decompositions [21 HD] • Now let A^, fc S N such that 
N /k S N, and assume that the MPS is injective (TU] such 
that for all j = fc, . . . , A — 2fc the set 



s- = 1 d 



(17) 



spans 
of 



5.,+fc+i ^ xhen \tp) is the unique ground state 



N/k-2 

71=0 



(18) 



where P„ is the projector onto the subspace orthogonal 
to the range of the mapping r„ : <cD„k+2k+iy<D„k+i _^ 



r„(A) = ^tr[AM„fe+i[ 



+ l] • • • Mnk+2k[Snk+ 



2k\ 



X \Snk+l ■ ■ ■ Snk+2k)- 

(19) 
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To get an efficiently computable lower bound the energy 
gap, we use 



AE = max <^ A H{H - A) > 



(20) 



We find 



p p I p p 



n,7n 



n,m 
l<|n-m|/c<2fc 



(21) 

where we omitted non-negative summands. Now consider 
the following quantity, which will bound the individual 
terms in the previous equation. 



7n,m min <{ A 
mm < A 



PnPm + PmPn + -^(^ri + Pm) ^ 
{Pn+Prr^f > (1- A)(P„+ An) 



= 1 — max < A 



(Pn + Pm)^ > XiPn+Pm) 



(22) 



where the maximum is given by the smallest non-zero 
eigenvalue of P„ + P,„ . Hence 

H^>H- J2 ln,mPn>{l-l)H, (23) 



and therefore we have the lower bound AE > 1 — 7, 
where 



7 — max 

0<n<iV 



E 



l<\7i-in\k<2k 



(24) 



Injectivity may fail to hold in certain singular cases. 
The simplest example given by the family of GHZ-type 
states ^(|0, . . . , 0) -f e''*|l, . . . , 1)). Since any local re- 
duced density matrix is independent of 0, it is impossible 
to distinguish the members of that family based on local 
information alone. Indeed, in this example, the ground 
state space of H will be two-dimensional, spanned by 
|0, . . . , 0) and |1, . . . , 1). One may check that the gap 
analysis above continues to hold in the degenerate case 
(unlike the original in [5]), now certifying the overlap 
between g and the ground-state space. This fact alone 
implies an exponential reduction in the number of un- 
known parameters. It is easy to see that the small re- 
maining uncertainty about g can be lifted in our example 
by measuring the "string operator" Gx ® ■ ■ ■ ® <Jx, which 
has expectation value cos Using e.g. the results of [38], 
the method of this particular example immediately gen- 
eralizes to any MPS with a two-fold degeneracy, such 
as the W state. Higher degeneracies may be treated by 
straight-forward, but more tedious, methods. 



l<\n-m\k<2k 
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